Run AlphaDiversity in scratchnotebooks
#source(here::here("RScripts", "InitialProcessing_3.R"))
library(vegan)
Loading required package: permute
Loading required package: lattice
This is vegan 2.6-3
library(cowplot)
library(flextable)
library(ftExtra)
source(here::here("RLibraries", "ChesapeakePersonalLibrary.R"))
This file is dedicated to conventional, non div-net/breakaway stats, since breakaway seems to choke on this data.
Reshape back into an ASV matrix, but this time correcting for total abundance
nonSpikes %>% head
raDf <- nonSpikes %>% pivot_wider(id_cols = ID, names_from = ASV, values_from = RA, values_fill = 0)
raMat <- raDf %>% column_to_rownames("ID")
raMat1 <- raMat %>% as.matrix()
countMat <- nonSpikes %>%
pivot_wider(id_cols = ID, names_from = ASV, values_from = reads, values_fill = 0) %>%
column_to_rownames("ID") %>% as.matrix()
seqDep <- countMat %>% apply(1, sum)
min(seqDep)
[1] 852
sampleRichness <- rarefy(countMat, min(seqDep))
rarefy everything to the minimum depth (852)
countRare <- rrarefy(countMat, min(seqDep))
Gamma diversity
specpool(countRare)
Doesn’t finish
#specpool(countMat)
All richness estimates
richnessRare <- estimateR(countRare)
Shannon diversity
shan <- diversity(countRare)
shan
3-1-B-0-2 3-1-B-1-2 3-1-B-180 3-1-B-20 3-1-B-5 3-1-B-500 3-1-B-53 3-1-S-0-2 3-1-S-1-2 3-1-S-180 3-1-S-20 3-1-S-5
4.444216 5.154575 4.645763 5.915606 5.132178 3.869332 5.560001 4.583410 4.902727 4.813431 5.336298 4.818610
3-2-B-0-2 3-2-B-1-2 3-2-B-180 3-2-B-20 3-2-B-5 3-2-B-500 3-2-B-53 3-2-S-0-2 3-2-S-1-2 3-2-S-180 3-2-S-20 3-2-S-5
4.344849 4.712513 4.675627 5.098950 5.442812 5.030542 4.890340 3.794814 4.924798 4.675124 4.995207 5.202696
3-2-S-500 3-2-S-53 3-3-B-0-2 3-3-B-1-2 3-3-B-180 3-3-B-20 3-3-B-5 3-3-B-500 3-3-B-53 3-3-S-180 3-3-S-20 3-3-S-500
4.883313 4.402086 4.437311 4.927824 3.459281 5.690676 5.258778 5.316690 5.477329 5.142843 4.917158 4.903164
3-3-S-53 4-3-B-0-2 4-3-B-1-2 4-3-B-180 4-3-B-20 4-3-B-5 4-3-B-500 4-3-B-53 4-3-O-1-2 4-3-O-180 4-3-O-5 4-3-O-500
4.361907 4.366410 4.965965 4.564515 4.452627 4.567276 4.231160 4.336562 5.019498 4.691336 5.147143 3.990599
4-3-O-53 4-3-S-0-2 4-3-S-180 4-3-S-20 4-3-S-500 4-3-S-53 5-1-S-1-2 5-1-S-180 5-1-S-20 5-1-S-5 5-1-S-500 5-1-S-53
4.550053 2.805419 4.479398 4.760915 4.733578 4.575868 4.459307 4.512310 4.168196 4.078014 4.011818 4.117664
5-5-B-0-2 5-5-B-180 5-5-B-5 5-5-B-500 5-5-B-53 5-5-S-180 5-5-S-5 5-5-S-500 5-5-S-53 C_5P1B_0P2 C_5P1B_180 C_5P1B_1P2
4.655036 5.151959 5.481992 5.093104 4.901198 4.294930 4.967999 4.952305 4.203941 4.368584 4.829209 4.897388
C_5P1B_20 C_5P1B_500 C_5P1B_53
5.468133 4.727174 5.189507
Evenness
pielouJ <- shan/richnessRare["S.chao1",]
pielouJ
3-1-B-0-2 3-1-B-1-2 3-1-B-180 3-1-B-20 3-1-B-5 3-1-B-500 3-1-B-53 3-1-S-0-2 3-1-S-1-2 3-1-S-180 3-1-S-20
0.010656840 0.006422546 0.011469251 0.002865951 0.006273652 0.041830621 0.007109976 0.008699531 0.007976329 0.008690465 0.007414011
3-1-S-5 3-2-B-0-2 3-2-B-1-2 3-2-B-180 3-2-B-20 3-2-B-5 3-2-B-500 3-2-B-53 3-2-S-0-2 3-2-S-1-2 3-2-S-180
0.009319248 0.011267761 0.009263164 0.012319590 0.006327012 0.007362969 0.008080497 0.010290423 0.020464414 0.008360626 0.010250272
3-2-S-20 3-2-S-5 3-2-S-500 3-2-S-53 3-3-B-0-2 3-3-B-1-2 3-3-B-180 3-3-B-20 3-3-B-5 3-3-B-500 3-3-B-53
0.007953267 0.008159670 0.008997352 0.015550953 0.012446876 0.007485125 0.061772871 0.005399056 0.008242407 0.006453177 0.007146198
3-3-S-180 3-3-S-20 3-3-S-500 3-3-S-53 4-3-B-0-2 4-3-B-1-2 4-3-B-180 4-3-B-20 4-3-B-5 4-3-B-500 4-3-B-53
0.008102995 0.010257642 0.010067304 0.013101847 0.010688464 0.009291003 0.009253490 0.009804310 0.008549470 0.010463207 0.009012522
4-3-O-1-2 4-3-O-180 4-3-O-5 4-3-O-500 4-3-O-53 4-3-S-0-2 4-3-S-180 4-3-S-20 4-3-S-500 4-3-S-53 5-1-S-1-2
0.009845773 0.012908904 0.008116955 0.012923854 0.008953551 0.124685269 0.011208241 0.008635149 0.009226702 0.012415445 0.009630661
5-1-S-180 5-1-S-20 5-1-S-5 5-1-S-500 5-1-S-53 5-5-B-0-2 5-5-B-180 5-5-B-5 5-5-B-500 5-5-B-53 5-5-S-180
0.012650427 0.013618132 0.015120557 0.012654152 0.012996712 0.008325702 0.008868699 0.007205562 0.009797355 0.009640319 0.013719828
5-5-S-5 5-5-S-500 5-5-S-53 C_5P1B_0P2 C_5P1B_180 C_5P1B_1P2 C_5P1B_20 C_5P1B_500 C_5P1B_53
0.009890166 0.010941140 0.016917268 0.007398312 0.008438195 0.005447460 0.003736643 0.006518782 0.005323390
diversityData <- sampleData %>%
left_join(richnessRare %>% t() %>% as.data.frame() %>% rownames_to_column("ID"), by = "ID") %>%
left_join(shan %>% enframe(name = "ID", value = "shannonH"), by = "ID") %>%
left_join(pielouJ %>% enframe(name = "ID", value = "pielouJ"), by = "ID") %>%
arrange(Size_Class)
Parameters for all plots
plotSpecs <- list(
facet_wrap(~Depth, ncol = 1) ,
theme_bw(base_size = 16) ,
geom_point(size = 4) ,
geom_path(aes(color = as.factor(Station))) ,
scale_x_log10(breaks = my_sizes, labels = as.character(my_sizes)) ,
#scale_y_log10nice() ,
scale_shape_manual(values = rep(21:25, 2)) ,
scale_fill_viridis_d(option = "plasma") ,
scale_color_viridis_d(option = "plasma") ,
labs(x = expression(paste("Particle Size (", mu, "m)"))) ,
theme(legend.position = "none",
plot.margin = unit(c(0, 0, 0, 0), "cm"),
axis.title.x = element_blank(),
axis.text.x = element_text(angle = 90, vjust = .5),
axis.title.y = element_text(margin = unit(c(3, 3, 3, 3), "mm"), vjust = 0))
)
Observed species counts, on rarefied data
plotObs <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = S.obs, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +ylab("Observed ASVs (Rarefied)")#+ scale_y_log10()
plotObs
Estemated Chao1 Richness
plotChao1 <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = S.chao1, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +
geom_errorbar(aes(ymin = S.chao1 -2 * se.chao1, ymax = S.chao1 + 2* se.chao1), width = -.1) +
scale_y_log10() +
ylab("Richness (Chao1)")
plotChao1
Shannon diversity
plotShan <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = shannonH, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +
ylab("Diversity (Shannon H)") +
lims(y = c(2.5, 6))
plotShan
Evenness
plotPielou <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pielouJ, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +scale_y_log10() +ylab("Evenness (PielouJ)")
plotPielou
All plots together
plotAlpha <- plot_grid(plotObs, plotChao1, plotShan, plotPielou, nrow = 1, labels = LETTERS)
plotAlpha
ggsave(here::here("Figures", "ConventionalAlpha.png"), plotAlpha, width = 11, height = 4)
Rarefied
obsMod <- lm(S.obs ~ log(Size_Class) + I(log(Size_Class)^2) + I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
summary(obsMod)
Call:
lm(formula = S.obs ~ log(Size_Class) + I(log(Size_Class)^2) +
I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
Residuals:
Min 1Q Median 3Q Max
-226.81 -37.14 4.33 44.82 206.47
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 283088.268 134955.453 2.098 0.039602 *
log(Size_Class) 30.953 8.177 3.786 0.000324 ***
I(log(Size_Class)^2) -6.124 1.521 -4.026 0.000144 ***
lat -14744.824 7032.354 -2.097 0.039687 *
I(lat^2) 192.061 91.548 2.098 0.039576 *
depth 4.929 3.238 1.522 0.132610
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 76.87 on 69 degrees of freedom
Multiple R-squared: 0.2527, Adjusted R-squared: 0.1986
F-statistic: 4.667 on 5 and 69 DF, p-value: 0.0009955
chao1Mod <- lm(S.chao1 ~ log(Size_Class) + I(log(Size_Class)^2) + I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
summary(chao1Mod)
Call:
lm(formula = S.chao1 ~ log(Size_Class) + I(log(Size_Class)^2) +
I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
Residuals:
Min 1Q Median 3Q Max
-546.23 -122.72 1.00 96.67 1310.85
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.031e+06 4.685e+05 2.202 0.03105 *
log(Size_Class) 7.652e+01 2.838e+01 2.696 0.00881 **
I(log(Size_Class)^2) -1.637e+01 5.281e+00 -3.101 0.00279 **
lat -5.378e+04 2.441e+04 -2.203 0.03095 *
I(lat^2) 7.008e+02 3.178e+02 2.205 0.03076 *
depth 2.267e+01 1.124e+01 2.016 0.04765 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 266.8 on 69 degrees of freedom
Multiple R-squared: 0.1903, Adjusted R-squared: 0.1316
F-statistic: 3.244 on 5 and 69 DF, p-value: 0.01091
chao1ModSimple <- lm(S.chao1 ~ log(Size_Class) + I(log(Size_Class)^2), data = diversityData)
summary(chao1ModSimple)
Call:
lm(formula = S.chao1 ~ log(Size_Class) + I(log(Size_Class)^2),
data = diversityData)
Residuals:
Min 1Q Median 3Q Max
-464.71 -138.96 -20.87 111.98 1409.67
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 574.310 48.667 11.801 < 2e-16 ***
log(Size_Class) 77.285 28.902 2.674 0.00927 **
I(log(Size_Class)^2) -16.870 5.376 -3.138 0.00247 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 272 on 72 degrees of freedom
Multiple R-squared: 0.1218, Adjusted R-squared: 0.09738
F-statistic: 4.992 on 2 and 72 DF, p-value: 0.009327
shanMod <- lm(shannonH ~ log(Size_Class) + I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
summary(shanMod)
Call:
lm(formula = shannonH ~ log(Size_Class) + I(log(Size_Class)^2) +
lat + I(lat^2) + depth, data = diversityData)
Residuals:
Min 1Q Median 3Q Max
-1.3356 -0.1482 0.0374 0.3250 0.7495
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.623e+03 7.861e+02 2.065 0.0427 *
log(Size_Class) 1.998e-01 4.763e-02 4.195 7.96e-05 ***
I(log(Size_Class)^2) -3.705e-02 8.861e-03 -4.182 8.36e-05 ***
lat -8.432e+01 4.096e+01 -2.059 0.0433 *
I(lat^2) 1.098e+00 5.332e-01 2.059 0.0433 *
depth 1.897e-02 1.886e-02 1.006 0.3180
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.4477 on 69 degrees of freedom
Multiple R-squared: 0.2971, Adjusted R-squared: 0.2461
F-statistic: 5.832 on 5 and 69 DF, p-value: 0.0001503
pielouMod <- lm(pielouJ ~ log(Size_Class) + I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityData)
summary(pielouMod)
Call:
lm(formula = pielouJ ~ log(Size_Class) + I(log(Size_Class)^2) +
lat + I(lat^2) + depth, data = diversityData)
Residuals:
Min 1Q Median 3Q Max
-0.014133 -0.004666 -0.002454 0.000611 0.100308
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -2.507e+01 2.635e+01 -0.952 0.3447
log(Size_Class) -3.775e-03 1.596e-03 -2.365 0.0209 *
I(log(Size_Class)^2) 6.540e-04 2.970e-04 2.202 0.0310 *
lat 1.306e+00 1.373e+00 0.951 0.3449
I(lat^2) -1.699e-02 1.787e-02 -0.950 0.3452
depth -3.932e-04 6.322e-04 -0.622 0.5361
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.01501 on 69 degrees of freedom
Multiple R-squared: 0.09128, Adjusted R-squared: 0.02544
F-statistic: 1.386 on 5 and 69 DF, p-value: 0.2402
uomisto H (2010a). “A diversity of beta diver- sities: straightening up a concept gone awry. 1. Defining beta diversity as a function of alpha and gamma diversity.” Ecography, 33, 2–2
predict(obsMod, se.fit = TRUE)
$fit
1 2 3 4 5 6 7 8 9 10 11 12 13 14
235.0858 235.0858 202.6115 202.6115 216.5492 165.3419 165.3419 196.1950 206.0323 306.2069 306.2069 273.7326 273.7326 287.6703
15 16 17 18 19 20 21 22 23 24 25 26 27 28
236.4630 236.4630 267.3161 267.3161 334.7205 334.7205 302.2462 302.2462 316.1839 264.9766 264.9766 295.8297 305.6670 305.6670
29 30 31 32 33 34 35 36 37 38 39 40 41 42
338.5321 338.5321 306.0579 306.0579 319.9956 319.9956 268.7882 268.7882 299.6413 299.6413 327.1209 294.6466 294.6466 308.5843
43 44 45 46 47 48 49 50 51 52 53 54 55 56
308.5843 257.3769 257.3769 257.3769 288.2300 288.2300 298.0673 298.0673 296.3515 296.3515 263.8772 263.8772 277.8150 277.8150
57 58 59 60 61 62 63 64 65 66 67 68 69 70
226.6076 226.6076 226.6076 257.4607 257.4607 267.2980 267.2980 256.5983 224.1240 224.1240 238.0617 238.0617 186.8544 186.8544
71 72 73 74 75
186.8544 217.7075 217.7075 227.5448 227.5448
$se.fit
[1] 26.89171 26.89171 26.15853 26.15853 26.01893 27.68111 27.68111 29.23541 33.91547 18.86417 18.86417 18.25002 18.25002 17.13889
[15] 19.96813 19.96813 21.42573 21.42573 19.18476 19.18476 18.72459 18.72459 17.17411 20.01466 20.01466 21.12886 28.07566 28.07566
[29] 19.18259 19.18259 18.69162 18.69162 16.90082 16.90082 19.54040 19.54040 20.50140 20.50140 18.49526 17.85851 17.85851 15.94069
[43] 15.94069 18.37681 18.37681 18.37681 19.36426 19.36426 26.45911 26.45911 19.24276 19.24276 18.35897 18.35897 16.62051 16.62051
[57] 18.36048 18.36048 18.36048 19.42883 19.42883 26.02108 26.02108 24.37188 23.42219 23.42219 22.26479 22.26479 23.05780 23.05780
[71] 23.05780 24.05365 24.05365 29.12182 29.12182
$df
[1] 69
$residual.scale
[1] 76.86895
diversityData$pred_obs = predict(obsMod, se.fit = TRUE)$fit
diversityData$se_obs = predict(obsMod, se.fit = TRUE)$se.fit
plotSpecs2 <- list(
facet_wrap(~Depth, ncol = 1) ,
theme_bw(base_size = 16) ,
#geom_point(size = 4) ,
geom_path(aes(color = as.factor(Station))) ,
scale_x_log10(breaks = my_sizes, labels = as.character(my_sizes)) ,
#scale_y_log10nice() ,
scale_shape_manual(values = rep(21:25, 2)) ,
scale_fill_viridis_d(option = "plasma") ,
scale_color_viridis_d(option = "plasma") ,
labs(x = expression(paste("Particle Size (", mu, "m)"))) ,
theme(legend.position = "none",
plot.margin = unit(c(0, 0, 0, 0), "cm"),
axis.title.x = element_blank(),
axis.text.x = element_text(angle = 90, vjust = .5),
axis.title.y = element_text(margin = unit(c(3, 3, 3, 3), "mm"), vjust = 0))
)
plotObs_pred <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pred_obs, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_obs - 2 * se_obs, yend = pred_obs + 2 * se_obs, xend = Size_Class, color = as.factor(Station)), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both"), alpha = 0.5) +
plotSpecs2 + ylab("Predicted ASVs")
plotObs_pred
predict(chao1Mod, se.fit = TRUE)
$fit
1 2 3 4 5 6 7 8 9 10 11 12 13 14
505.4108 505.4108 379.2309 379.2309 480.8810 295.0839 295.0839 434.4257 374.1099 684.3823 684.3823 558.2024 558.2024 659.8525
15 16 17 18 19 20 21 22 23 24 25 26 27 28
474.0554 474.0554 613.3972 613.3972 751.7117 751.7117 625.5318 625.5318 727.1819 541.3848 541.3848 680.7266 620.4108 620.4108
29 30 31 32 33 34 35 36 37 38 39 40 41 42
753.2514 753.2514 627.0715 627.0715 728.7216 728.7216 542.9245 542.9245 682.2663 682.2663 716.6596 590.4797 590.4797 692.1298
43 44 45 46 47 48 49 50 51 52 53 54 55 56
692.1298 506.3326 506.3326 506.3326 645.6744 645.6744 585.3586 585.3586 626.7627 626.7627 500.5828 500.5828 602.2329 602.2329
57 58 59 60 61 62 63 64 65 66 67 68 69 70
416.4358 416.4358 416.4358 555.7776 555.7776 495.4618 495.4618 514.1003 387.9204 387.9204 489.5705 489.5705 303.7734 303.7734
71 72 73 74 75
303.7734 443.1152 443.1152 382.7994 382.7994
$se.fit
[1] 93.34957 93.34957 90.80446 90.80446 90.31987 96.08980 96.08980 101.48527 117.73123 65.48346 65.48346 63.35153
[13] 63.35153 59.49448 69.31563 69.31563 74.37542 74.37542 66.59631 66.59631 64.99893 64.99893 59.61672 69.47716
[25] 69.47716 73.34491 97.45941 97.45941 66.58878 66.58878 64.88447 64.88447 58.66804 58.66804 67.83086 67.83086
[37] 71.16680 71.16680 64.20284 61.99249 61.99249 55.33512 55.33512 63.79168 63.79168 63.79168 67.21942 67.21942
[49] 91.84787 91.84787 66.79767 66.79767 63.72976 63.72976 57.69500 57.69500 63.73497 63.73497 63.73497 67.44355
[61] 67.44355 90.32734 90.32734 84.60242 81.30576 81.30576 77.28806 77.28806 80.04086 80.04086 80.04086 83.49774
[73] 83.49774 101.09096 101.09096
$df
[1] 69
$residual.scale
[1] 266.8362
diversityData$pred_chao1 = predict(chao1Mod, se.fit = TRUE)$fit
diversityData$se_chao1 = predict(chao1Mod, se.fit = TRUE)$se.fit
plotChao1_pred <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pred_chao1, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_chao1 - 2 * se_chao1, yend = pred_chao1 + 2 * se_chao1, xend = Size_Class, color = as.factor(Station)), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both"), alpha = 0.5) +
plotSpecs2 + ylab("Predictd Richness (Chao1)") + scale_y_log10()
plotChao1_pred
predict(shanMod, se.fit = TRUE)
$fit
1 2 3 4 5 6 7 8 9 10 11 12 13 14
4.531951 4.531951 4.360594 4.360594 4.338910 4.039890 4.039890 4.161410 4.400479 4.984722 4.984722 4.813365 4.813365 4.791681
15 16 17 18 19 20 21 22 23 24 25 26 27 28
4.492661 4.492661 4.614181 4.614181 5.175134 5.175134 5.003777 5.003777 4.982092 4.683073 4.683073 4.804592 5.043662 5.043662
29 30 31 32 33 34 35 36 37 38 39 40 41 42
5.215580 5.215580 5.044223 5.044223 5.022539 5.022539 4.723519 4.723519 4.845039 4.845039 5.158760 4.987403 4.987403 4.965719
43 44 45 46 47 48 49 50 51 52 53 54 55 56
4.965719 4.666699 4.666699 4.666699 4.788219 4.788219 5.027288 5.027288 4.987930 4.987930 4.816573 4.816573 4.794889 4.794889
57 58 59 60 61 62 63 64 65 66 67 68 69 70
4.495869 4.495869 4.495869 4.617389 4.617389 4.856458 4.856458 4.760224 4.588867 4.588867 4.567183 4.567183 4.268163 4.268163
71 72 73 74 75
4.268163 4.389683 4.389683 4.628752 4.628752
$se.fit
[1] 0.15663934 0.15663934 0.15236870 0.15236870 0.15155556 0.16123742 0.16123742 0.17029094 0.19755145 0.10988038 0.10988038
[12] 0.10630304 0.10630304 0.09983095 0.11631071 0.11631071 0.12480098 0.12480098 0.11174773 0.11174773 0.10906736 0.10906736
[23] 0.10003607 0.11658176 0.11658176 0.12307179 0.16353561 0.16353561 0.11173510 0.11173510 0.10887529 0.10887529 0.09844420
[34] 0.09844420 0.11381929 0.11381929 0.11941696 0.11941696 0.10773152 0.10402259 0.10402259 0.09285161 0.09285161 0.10704160
[45] 0.10704160 0.10704160 0.11279330 0.11279330 0.15411952 0.15411952 0.11208560 0.11208560 0.10693769 0.10693769 0.09681145
[56] 0.09681145 0.10694644 0.10694644 0.10694644 0.11316940 0.11316940 0.15156809 0.15156809 0.14196175 0.13643000 0.13643000
[67] 0.12968835 0.12968835 0.13430751 0.13430751 0.13430751 0.14010812 0.14010812 0.16962930 0.16962930
$df
[1] 69
$residual.scale
[1] 0.4477477
diversityData$pred_shanH = predict(shanMod, se.fit = TRUE)$fit
diversityData$se_shanH = predict(shanMod, se.fit = TRUE)$se.fit
plotShannonH_pred <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pred_shanH, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_shanH - 2 * se_shanH, yend = pred_shanH + 2 * se_shanH, xend = Size_Class, color = as.factor(Station)), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both"), alpha = 0.5) +
plotSpecs2 + ylab("Predicted Diversity (Shannon H)") #+ scale_y_log10()
plotShannonH_pred
predict(pielouMod, se.fit = TRUE)
$fit
1 2 3 4 5 6 7 8 9 10 11
0.018866864 0.018866864 0.021537876 0.021537876 0.020512033 0.024377201 0.024377201 0.021531485 0.018827124 0.010430513 0.010430513
12 13 14 15 16 17 18 19 20 21 22
0.013101525 0.013101525 0.012075681 0.015940850 0.015940850 0.013095134 0.013095134 0.006715502 0.006715502 0.009386514 0.009386514
23 24 25 26 27 28 29 30 31 32 33
0.008360671 0.012225839 0.012225839 0.009380123 0.006675762 0.006675762 0.005657637 0.005657637 0.008328649 0.008328649 0.007302806
34 35 36 37 38 39 40 41 42 43 44
0.007302806 0.011167974 0.011167974 0.008322258 0.008322258 0.006418762 0.009089774 0.009089774 0.008063930 0.008063930 0.011929098
45 46 47 48 49 50 51 52 53 54 55
0.011929098 0.011929098 0.009083382 0.009083382 0.006379021 0.006379021 0.009130701 0.009130701 0.011801713 0.011801713 0.010775870
56 57 58 59 60 61 62 63 64 65 66
0.010775870 0.014641038 0.014641038 0.014641038 0.011795322 0.011795322 0.009090960 0.009090960 0.012896424 0.015567436 0.015567436
67 68 69 70 71 72 73 74 75
0.014541593 0.014541593 0.018406761 0.018406761 0.018406761 0.015561045 0.015561045 0.012856683 0.012856683
$se.fit
[1] 0.005250060 0.005250060 0.005106922 0.005106922 0.005079668 0.005404174 0.005404174 0.005707619 0.006621306 0.003682846
[11] 0.003682846 0.003562945 0.003562945 0.003346021 0.003898371 0.003898371 0.004182938 0.004182938 0.003745434 0.003745434
[21] 0.003655596 0.003655596 0.003352896 0.003907456 0.003907456 0.004124981 0.005481202 0.005481202 0.003745011 0.003745011
[31] 0.003649159 0.003649159 0.003299541 0.003299541 0.003814866 0.003814866 0.004002483 0.004002483 0.003610823 0.003486512
[41] 0.003486512 0.003112095 0.003112095 0.003587699 0.003587699 0.003587699 0.003780478 0.003780478 0.005165604 0.005165604
[51] 0.003756759 0.003756759 0.003584217 0.003584217 0.003244817 0.003244817 0.003584510 0.003584510 0.003584510 0.003793084
[61] 0.003793084 0.005080088 0.005080088 0.004758113 0.004572707 0.004572707 0.004346748 0.004346748 0.004501567 0.004501567
[71] 0.004501567 0.004695985 0.004695985 0.005685443 0.005685443
$df
[1] 69
$residual.scale
[1] 0.0150071
diversityData$pred_pielouJ = predict(pielouMod, se.fit = TRUE)$fit
diversityData$se_pielouJ = predict(pielouMod, se.fit = TRUE)$se.fit
plot_pielouJ_pred <- diversityData %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pred_pielouJ, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_pielouJ - 2 * se_pielouJ, yend = pred_pielouJ + 2 * se_pielouJ, xend = Size_Class, color = as.factor(Station), alpha = 0.5), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both")) +
plotSpecs2 + ylab("Predicted Evenness (Pielou J)") + scale_y_log10()
plot_pielouJ_pred
plotPredictions <- plot_grid(plotObs_pred, plotChao1_pred, plotShannonH_pred, plot_pielouJ_pred, nrow = 1, labels = LETTERS)
Warning: NaNs producedWarning: Transformation introduced infinite values in continuous y-axisWarning: Removed 11 rows containing missing values (geom_segment).
plotPredictions
ggsave(here::here("Figures", "ConventionalAlphaPredictions.png"), plotPredictions, width = 11, height = 4)
plot_grid(plotObs, plotChao1, plotShan, plotPielou,
plotObs_pred, plotChao1_pred, plotShannonH_pred, plot_pielouJ_pred,
nrow = 2, labels = LETTERS)
Warning: NaNs producedWarning: Transformation introduced infinite values in continuous y-axisWarning: Removed 11 rows containing missing values (geom_segment).
alphaSummary <- tibble(
metric = c("Observed ASVs", "Richness (Chao1)", "Diversity (Shannon H)", "Evenness (Pielou J)"),
model = list(obsMod, chao1Mod, shanMod, pielouMod)
)
alphaSummary <- alphaSummary %>%
mutate(df = map(model, ~broom::tidy(summary(.))))
alphaSummary <- alphaSummary %>%
select(-model) %>%
unnest(df)
alphaSummary <- alphaSummary %>%
rename(Metric = metric, Term = term, Estimate = estimate, `Standard Error` = std.error, `T Value` = statistic, p = p.value) %>%
mutate(Term = str_replace(Term, "^I?\\((.*)\\)", "\\1"),
Term = str_replace(Term, "\\^2", "\\^2\\^"),
Term = str_replace(Term, "depth", "Depth"),
Term = str_replace(Term, "lat", "Latitude"),
Term = str_replace(Term, "_", " ")# BOOKMARK!!
) %>%
mutate(Estimate = format(Estimate, digits = 2, scientific = TRUE) %>%
reformat_sci()
) %>%
mutate(`Standard Error` = format(`Standard Error`, digits = 2, scientific = TRUE) %>%
reformat_sci()
) %>%
mutate(`T Value` = format(`T Value`, digits = 2, scientific = FALSE)) %>%
mutate(p = if_else(p < 0.001, "< 0.001", format(round(p, digits = 3)))) %>%
rename(`Standard\nError` = `Standard Error`) %>%
identity()
alphaSummary %>% flextable() %>% merge_v(j = 1) %>% theme_vanilla() %>% bold(i = ~ p< 0.05, j = "p") %>% colformat_md() %>% set_table_properties(layout = "autofit") %>%
align(j = -c(1:2), align = "right")
Metric | Term | Estimate | Standard | T Value | p |
Observed ASVs | Intercept | 2.8 x 105 | 1.3 x 105 | 2.10 | 0.040 |
log(Size Class) | 3.1 x 101 | 8.2 x 100 | 3.79 | < 0.001 | |
log(Size Class)2 | -6.1 x 100 | 1.5 x 100 | -4.03 | < 0.001 | |
Latitude | -1.5 x 104 | 7.0 x 103 | -2.10 | 0.040 | |
Latitude2 | 1.9 x 102 | 9.2 x 101 | 2.10 | 0.040 | |
Depth | 4.9 x 100 | 3.2 x 100 | 1.52 | 0.133 | |
Richness (Chao1) | Intercept | 1.0 x 106 | 4.7 x 105 | 2.20 | 0.031 |
log(Size Class) | 7.7 x 101 | 2.8 x 101 | 2.70 | 0.009 | |
log(Size Class)2 | -1.6 x 101 | 5.3 x 100 | -3.10 | 0.003 | |
Latitude | -5.4 x 104 | 2.4 x 104 | -2.20 | 0.031 | |
Latitude2 | 7.0 x 102 | 3.2 x 102 | 2.21 | 0.031 | |
Depth | 2.3 x 101 | 1.1 x 101 | 2.02 | 0.048 | |
Diversity (Shannon H) | Intercept | 1.6 x 103 | 7.9 x 102 | 2.06 | 0.043 |
log(Size Class) | 2.0 x 10-1 | 4.8 x 10-2 | 4.20 | < 0.001 | |
log(Size Class)2 | -3.7 x 10-2 | 8.9 x 10-3 | -4.18 | < 0.001 | |
Latitude | -8.4 x 101 | 4.1 x 101 | -2.06 | 0.043 | |
Latitude2 | 1.1 x 100 | 5.3 x 10-1 | 2.06 | 0.043 | |
Depth | 1.9 x 10-2 | 1.9 x 10-2 | 1.01 | 0.318 | |
Evenness (Pielou J) | Intercept | -2.5 x 101 | 2.6 x 101 | -0.95 | 0.345 |
log(Size Class) | -3.8 x 10-3 | 1.6 x 10-3 | -2.36 | 0.021 | |
log(Size Class)2 | 6.5 x 10-4 | 3.0 x 10-4 | 2.20 | 0.031 | |
Latitude | 1.3 x 100 | 1.4 x 100 | 0.95 | 0.345 | |
Latitude2 | -1.7 x 10-2 | 1.8 x 10-2 | -0.95 | 0.345 | |
Depth | -3.9 x 10-4 | 6.3 x 10-4 | -0.62 | 0.536 |
richSummary %>% rename_(.dots = setNames(names(.), paste0('break_', names(.))))
Warning: `rename_()` was deprecated in dplyr 0.7.0.
Please use `rename()` instead.
diversityDataWB <- full_join(diversityData,
richSummary %>% rename_(.dots = setNames(names(.), paste0('break_', names(.)))),
by = c("ID" = "break_sample_names"), suffix = c("", "_break")) %>%
mutate(pielouJ2 = shannonH/break_estimate) %>%
identity()
diversityDataWB
pielouMod2 <- lm(pielouJ2 ~ log(Size_Class) + I(log(Size_Class)^2) + lat + I(lat^2) + depth, data = diversityDataWB)
summary(pielouMod2)
Call:
lm(formula = pielouJ2 ~ log(Size_Class) + I(log(Size_Class)^2) +
lat + I(lat^2) + depth, data = diversityDataWB)
Residuals:
Min 1Q Median 3Q Max
-0.013938 -0.005151 -0.002547 0.000861 0.105939
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.962e+01 2.755e+01 -0.712 0.4788
log(Size_Class) -3.288e-03 1.669e-03 -1.970 0.0529 .
I(log(Size_Class)^2) 5.742e-04 3.106e-04 1.849 0.0688 .
lat 1.020e+00 1.436e+00 0.711 0.4797
I(lat^2) -1.325e-02 1.869e-02 -0.709 0.4806
depth -2.376e-04 6.612e-04 -0.359 0.7204
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.01569 on 69 degrees of freedom
Multiple R-squared: 0.06714, Adjusted R-squared: -0.0004537
F-statistic: 0.9933 on 5 and 69 DF, p-value: 0.4284
Ok. So the narrative makes sense. Alpha diveristy is driven by variability in richness rather than evenness. Why would we see an effect in chao1 but not breakaway? Because chao1 is more driven by abundant stuff that makes the rarification threshold. My first hunch is that chao1 responds to evenness, but actually that shouldn’t have any effect since there is no evenness variability? Or maybe just that breakaway is more variable (because it detects fine level differences in rare species) and that doesn’t map as nicely with overall patterns.
plotBreak <- diversityDataWB %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = break_estimate, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +
#scale_y_log10()+
ylab("Richness (Breakaway)")
plotBreak
plotPielou2 <- diversityDataWB %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pielouJ2, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +
#scale_y_log10()+
ylab("Evenness (PielouJ)")
plotPielou2
predict(pielouMod2, se.fit = TRUE)
$fit
1 2 3 4 5 6 7 8 9
0.0117020065 0.0117020065 0.0135963432 0.0135963432 0.0133971262 0.0160361854 0.0160361854 0.0140005291 0.0099098120
10 11 12 13 14 15 16 17 18
0.0043418030 0.0043418030 0.0062361398 0.0062361398 0.0060369228 0.0086759820 0.0086759820 0.0066403257 0.0066403257
19 20 21 22 23 24 25 26 27
0.0011173144 0.0011173144 0.0030116511 0.0030116511 0.0028124341 0.0054514933 0.0054514933 0.0034158370 -0.0006748801
28 29 30 31 32 33 34 35 36
-0.0006748801 0.0002246683 0.0002246683 0.0021190051 0.0021190051 0.0019197880 0.0019197880 0.0045588472 0.0045588472
37 38 39 40 41 42 43 44 45
0.0025231909 0.0025231909 0.0009183182 0.0028126550 0.0028126550 0.0026134379 0.0026134379 0.0052524972 0.0052524972
46 47 48 49 50 51 52 53 54
0.0052524972 0.0032168409 0.0032168409 -0.0008738762 -0.0008738762 0.0033312013 0.0033312013 0.0052255380 0.0052255380
55 56 57 58 59 60 61 62 63
0.0050263210 0.0050263210 0.0076653802 0.0076653802 0.0076653802 0.0056297239 0.0056297239 0.0015390069 0.0015390069
64 65 66 67 68 69 70 71 72
0.0066640405 0.0085583773 0.0085583773 0.0083591603 0.0083591603 0.0109982195 0.0109982195 0.0109982195 0.0089625632
73 74 75
0.0089625632 0.0048718461 0.0048718461
$se.fit
[1] 0.005490363 0.005490363 0.005340673 0.005340673 0.005312171 0.005651530 0.005651530 0.005968865 0.006924373 0.003851415
[11] 0.003851415 0.003726026 0.003726026 0.003499173 0.004076805 0.004076805 0.004374397 0.004374397 0.003916868 0.003916868
[21] 0.003822918 0.003822918 0.003506363 0.004086305 0.004086305 0.004313787 0.005732084 0.005732084 0.003916425 0.003916425
[31] 0.003816186 0.003816186 0.003450566 0.003450566 0.003989478 0.003989478 0.004185682 0.004185682 0.003776096 0.003646094
[41] 0.003646094 0.003254540 0.003254540 0.003751913 0.003751913 0.003751913 0.003953516 0.003953516 0.005402041 0.005402041
[51] 0.003928711 0.003928711 0.003748271 0.003748271 0.003393337 0.003393337 0.003748578 0.003748578 0.003748578 0.003966699
[61] 0.003966699 0.005312611 0.005312611 0.004975899 0.004782006 0.004782006 0.004545704 0.004545704 0.004707610 0.004707610
[71] 0.004707610 0.004910927 0.004910927 0.005945674 0.005945674
$df
[1] 69
$residual.scale
[1] 0.015694
diversityDataWB$pred_pielouJ2 = predict(pielouMod2, se.fit = TRUE)$fit
diversityDataWB$se_pielouJ2 = predict(pielouMod2, se.fit = TRUE)$se.fit
plot_pielouJ2_pred <- diversityDataWB %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = pred_pielouJ2, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_pielouJ2 - 2 * se_pielouJ2, yend = pred_pielouJ2 + 2 * se_pielouJ2, xend = Size_Class, color = as.factor(Station), alpha = 0.5), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both")) +
plotSpecs2 + ylab("Predicted Evenness (Pielou J2)") #+ scale_y_log10()
plot_pielouJ2_pred
plotBreakaway <- diversityDataWB %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
ggplot(aes(x = Size_Class, y = break_estimate, shape = as.factor(Station), fill = as.factor(Station))) +
plotSpecs +
geom_errorbar(aes(ymin = break_lower, ymax = break_upper), width = -.1) +
scale_y_log10() +
ylab("Richness (Breakaway)")
plotBreakaway
#predict(breakLm, se.fit = TRUE)
# doesn't work because built with a different data frame
Why are these not smooth curves?!! What if I redo the model, this time with the same data frame
breakLm2 <- lm(break_estimate ~ log(Size_Class) + I(log(Size_Class) ^2) + lat + I(lat^2) + depth ,data = diversityDataWB)
breakLm2 %>% summary()
Call:
lm(formula = break_estimate ~ log(Size_Class) + I(log(Size_Class)^2) +
lat + I(lat^2) + depth, data = diversityDataWB)
Residuals:
Min 1Q Median 3Q Max
-2974.5 -1191.2 -151.6 599.9 6768.1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7124615.61 3339862.88 2.133 0.0365 *
log(Size_Class) 244.45 202.35 1.208 0.2312
I(log(Size_Class)^2) -75.16 37.65 -1.996 0.0498 *
lat -370568.38 174035.93 -2.129 0.0368 *
I(lat^2) 4817.28 2265.61 2.126 0.0371 *
depth 151.10 80.15 1.885 0.0636 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1902 on 69 degrees of freedom
Multiple R-squared: 0.1414, Adjusted R-squared: 0.0792
F-statistic: 2.273 on 5 and 69 DF, p-value: 0.0567
Note the non statistical significance overall
#predict(breakLm2, se.fit = TRUE)
diversityDataWB$pred_break = predict(breakLm2, se.fit = TRUE)$fit
diversityDataWB$se_break = predict(breakLm2, se.fit = TRUE)$se.fit
plot_break_pred <- diversityDataWB %>%
filter(Depth %in% c("Surface", "Bottom")) %>%
# filter(Station == 4.3) %>%
ggplot(aes(x = Size_Class, y = pred_break, shape = as.factor(Station), fill = as.factor(Station))) +
geom_segment(aes(y = pred_break - 2 * se_break, yend = pred_break + 2 * se_break, xend = Size_Class, color = as.factor(Station), alpha = 0.5), arrow = arrow(angle = 70, length = unit(0.05, "in"), ends = "both")) +
plotSpecs2 + ylab("Predicted Richness (Breakaway -- LM)") #+ scale_y_log10()
plot_break_pred
plotAlphaWB <- plot_grid(plotBreakaway, plotShan, plotPielou2, nrow = 1, labels = LETTERS)
plotAlphaWB
ggsave(here::here("Figures", "BreakawayAlpha.png"), plotAlpha, width = 11, height = 4)
Summary table I want both breakaway metrics here
bettaTable <- myBet$table %>%
as.data.frame() %>%
rename(estimate = Estimates,
`std.error` = `Standard Errors`,
`p.value`=`p-values`
) %>%
mutate(`statistic` = NA) %>%
rownames_to_column(var = "term") %>%
select(term, estimate, std.error, statistic, p.value) %>%
as_tibble()
bettaTable
alphaSummary2 <- tibble(
metric = c("Richness (Breakaway -- LM)", "Diversity (Shannon H)", "Evenness (Pielou J)"),
model = list(breakLm, shanMod, pielouMod2)
)
alphaSummary2 <- alphaSummary2 %>%
mutate(df = map(model, ~broom::tidy(summary(.))))
## Add in willis variables
breakawaySummary <- tibble(
metric = "Richness (Breakaway -- Betta)",
model = NULL,
df = list(bettaTable)
)
alphaSummary2 = bind_rows(breakawaySummary, alphaSummary2)
alphaSummary2 <- alphaSummary2 %>%
select(-model) %>%
unnest(df)
alphaSummary2 <- alphaSummary2 %>%
rename(Metric = metric, Term = term, Estimate = estimate, `Standard Error` = std.error, `T Value` = statistic, p = p.value) %>%
mutate(Term = str_replace(Term, "^I?\\((.*)\\)", "\\1"),
Term = str_replace(Term, "\\^2", "\\^2\\^"),
Term = str_replace(Term, "depth", "Depth"),
Term = str_replace(Term, "lat", "Latitude"),
Term = str_replace(Term, "_", " ")# BOOKMARK!!
) %>%
mutate(Estimate = format(Estimate, digits = 2, scientific = TRUE) %>%
reformat_sci()
) %>%
mutate(`Standard Error` = format(`Standard Error`, digits = 2, scientific = TRUE) %>%
reformat_sci()
) %>%
mutate(`T Value` = format(`T Value`, digits = 2, scientific = FALSE)) %>%
mutate(p = if_else(p < 0.001, "< 0.001", format(round(p, digits = 3)))) %>%
rename(`Standard\nError` = `Standard Error`) %>%
identity()
alphaSummary2
alphaTable2 <- alphaSummary2 %>% flextable() %>% merge_v(j = 1) %>% theme_vanilla() %>% bold(i = ~ p< 0.05, j = "p") %>% colformat_md() %>% set_table_properties(layout = "autofit") %>%
align(j = -c(1:2), align = "right")
alphaTable2
Metric | Term | Estimate | Standard | T Value | p |
Richness (Breakaway – Betta) | Intercept | 7.1 x 106 | 2.4 x 102 | NA | < 0.001 |
log(Size Class) | 1.2 x 102 | 6.1 x 101 | NA | 0.058 | |
log(Size Class)2 | -5.0 x 101 | 1.2 x 101 | NA | < 0.001 | |
Latitude | -3.7 x 105 | 6.1 x 100 | NA | < 0.001 | |
Latitude2 | 4.8 x 103 | 1.6 x 10-1 | NA | < 0.001 | |
Depth | 1.5 x 102 | 1.0 x 101 | NA | < 0.001 | |
Richness (Breakaway – LM) | Intercept | 7.1 x 106 | 3.3 x 106 | 2.13 | 0.036 |
log(Size Class) | 2.4 x 102 | 2.0 x 102 | 1.21 | 0.231 | |
log(Size Class)2 | -7.5 x 101 | 3.8 x 101 | -2.00 | 0.050 | |
Latitude | -3.7 x 105 | 1.7 x 105 | -2.13 | 0.037 | |
Latitude2 | 4.8 x 103 | 2.3 x 103 | 2.13 | 0.037 | |
Depth | 1.5 x 102 | 8.0 x 101 | 1.89 | 0.064 | |
Diversity (Shannon H) | Intercept | 1.6 x 103 | 7.9 x 102 | 2.06 | 0.043 |
log(Size Class) | 2.0 x 10-1 | 4.8 x 10-2 | 4.20 | < 0.001 | |
log(Size Class)2 | -3.7 x 10-2 | 8.9 x 10-3 | -4.18 | < 0.001 | |
Latitude | -8.4 x 101 | 4.1 x 101 | -2.06 | 0.043 | |
Latitude2 | 1.1 x 100 | 5.3 x 10-1 | 2.06 | 0.043 | |
Depth | 1.9 x 10-2 | 1.9 x 10-2 | 1.01 | 0.318 | |
Evenness (Pielou J) | Intercept | -2.0 x 101 | 2.8 x 101 | -0.71 | 0.479 |
log(Size Class) | -3.3 x 10-3 | 1.7 x 10-3 | -1.97 | 0.053 | |
log(Size Class)2 | 5.7 x 10-4 | 3.1 x 10-4 | 1.85 | 0.069 | |
Latitude | 1.0 x 100 | 1.4 x 100 | 0.71 | 0.480 | |
Latitude2 | -1.3 x 10-2 | 1.9 x 10-2 | -0.71 | 0.481 | |
Depth | -2.4 x 10-4 | 6.6 x 10-4 | -0.36 | 0.720 |
alphaTable2 %>% save_as_docx(path = here::here("Tables", "alphaTable2.docx"))
myBet$table
plot_grid(plot_break_pred,plotShannonH_pred,plot_pielouJ2_pred, nrow = 1, labels = LETTERS)